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In this paper, a memristive Murah-Lakshmanan-Chua (MLC) circuit is built by replacing the 
nonlinear element of an ordinary MLC circuit, namely the Chua's diode, with a three segment 
piecewise linear active flux controlled memristor. The bist ability nature of the memristor introduces 
two discontinuty boundaries or switching manifolds in the circuit topology. As a result, the circuit 
becomes a piecewise smooth system of second order. Grazing bifurcations, which are essentially a 
form of discontinuity induced non-smooth bifurcations, occur at these boundaries and govern the 
dynamics of the circuit. While the interaction of the memristor aided self oscillations of the circuit 
and the external sinusoidal forcing result in the phenomenon of beats occurring in the circuit, 
grazing bifurcations endow them with chaotic and hyper chaotic nature. In addition the circuit 
admits a codimension-5 bifurcation and transient hyper chaos. Grazing bifurcations as well as 
other behaviors have been analyzed numerically using time series plots, phase portraits, bifurcation 
diagram, power spectra and Lyapunov spectrum, as well as the recent 0-1 K test for chaos, obtained 
after constructing a proper Zero Time Discontinuity Map (ZDM) and Poincare Discontinuity Map 
(PDM) analytically. Multisim simulations using a model of piecewise linear memristor have also 
been used to confirm some of the behaviors. 



I. INTRODUCTION 

The Murali-Lakshmanan-Chua (MLC) circuit which 
is a two dimensional circuit first introduced by Murali 
et al. [1 , is basically a classic configuration of forced se- 
ries LCR oscillatory circuit having a Chua's diode as its 
non linearity. It is found to exhibit a large variety of 
bifurcation and chaos phenomena [T]-[3]. An eigen value 
study of this circuit was made by Lindberg and Murali 
[4 . An exhaustive experimental study on it has been per- 
formed by Lakshmanan and Murali [5 and on the MLC 
and a variant of MLC circuits by Thamilmaran and Lak- 
shmanan [HIT]- The control and synchronization of chaos 
in this circuit has been effected by different mechanisms 
such as nonfeedback methods [StiTT] and in the presence 
of noise by Ramesh and Narayanan [12]. The spatiotem- 
poral dynamics of coupled array of this circuit has been 
studied by Muruganandam et al. [13 . The birth of SNA 
through type III intermittency route has been reported 
by Venkatesan et al. [14 . While a statistical characteri- 
zation of the chaotic attractors at bifurcations has been 
carried out by Philominathan and Neelamegam [15 , a 
new class of chaotic attractors in this circuit has been re- 
ported by Oksasoglu and Wang [T6|. The observation of 



chaotic beats in a quasi periodically forced MLC circuit 
has been reported in [17 . An inductorless realization of 
this circuit has been designed using current feedback op- 
erational amplifiers (CFOA) by Cam and Kilic [T8^ and 
a mixed mode operation of this circuit using CFOA's has 
been realized by Cam [19] and Kilic [20]. The presence 
of multiple period doubling bifurcation route to chaos 
when periodically pulsed perturbations are applied has 
been identified by Srinivasan [21 . A theoretical study of 
the memristive MLC circuit has been reported by Wang 
et al. [22 . Chaos in fractional order MLC circuit has 
been studied in [23 and in a fractional memristor based 
MLC system with cubic non linearity by [24 . The MLC 
circuit in the frame of CNN has been studied by Gunay 
[25]. The ordered and chaotic dynamics of two coupled 
MLC like series LCR circuits with a common non lin- 
earity has been studied by Santhiah et al. [26 . Thus 
we find that over the years the individual, coupled and 
spatio-temporal dynamics of the MLC circuit have been 
intensively investigated by a large number of workers. Its 
importance lies in its low dimensionality, conceptual sim- 
plicity as well as the mathematically tractable nature of 
its nonlinear element - the Chua's diode. 
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A. Motivation and plan 

The motivation for the present work is the fohowing. 

• Is it possible to modify the MLC circuit into a mem- 
ristive MLC circuit by the addition / replacement 
of its non linearity with a multisim prototype model 
of the memristor introduced by the present authors 
[27] in an earlier work on Chua's circuit? 

• Will the addition of this memristor make the re- 
sultant circuit a non-smooth or piecewise-smooth 
dynamical system? 

• If so, what are the types of non- smooth bifurcations 
that the circuit may admit? 

• Will the addition of memristor in the MLC circuit 
result in the same type of switching action and oc- 
currence of beats that were observed in the mem- 
ristive Chua's circuit? 

• If so, what will be the nature of the beats: 
quasiperiodic or chaotic or hyper chaotic? 

The following work is very much directed towards find- 
ing answers to these questions. We have built a memris- 
tive MLC circuit, similar to the one reported by Wang 
et al. [22 numerically, by replacing the Chua's diode with 
a three segment piecewise fiux controlled active memris- 
tor deisgned by the same authors earlier, as its non lin- 
earity. We have found that the addition of memristor 
as the nonlinear element converts the system into a sec- 
ond order piece- wise continuous system having two dis- 
continuous boundaries which admit grazing bifurcations. 
The non-smooth dynamics and the time varying resistive 
property and the switching characteristic of the memris- 
tor give rise to many interesting behaviors such as graz- 
ing bifurcations, hyper chaos, transient hyperchaos and 
hyper chaotic beats in this circuit for a chosen sets of 
circuit parameters. To our knowledge, we believe that it 
is for the first time that hyper chaotic beats have been 
generated in a simple system, non-smooth or otherwise, 
having just a single driving force. 

The plan of the paper is as follows. In Sec. 2, the 
characteristic of a Chua-type fiux controlled active mem- 
ristor and its switching properties are discussed and a 
prototype Multisim model of it, developed by the same 
authors in an earlier work, is described. In Sec. 3, a 
modification of the standard MLC circuit by addition of 
this fiux controlled memristor is discussed and the cir- 
cuit equations are written based upon Kirchoff's laws. 
The stability analysis of this circuit is also made here. 
In Sec. 4 a general introduction to non-smooth systems 
and the types of bifurcations admitted by such systems 
and the means of observing them correctly are given. In 
Sec. 5 the classification of the memristive MLC circuit 
as a nonsmooth system of second order admitting graz- 
ing bifurcations is discussed. Further the construction 
of the Zero Discontinuity Map (ZDM) and the Poincare 



Discontinuity Map (PDM) are also mentioned. In Sec. 
6 numerical observation of grazing bifurcation is given. 
That grazing bifurcations lead to a change in the dynam- 
ics of the system has been verified through the applica- 
tion of the 0-1 binary or the K test for chaos, in addition 
to the usual verification using Lyapunov exponents. The 
codimension-5 bifurcation observed in this circuit is re- 
ported in Sec. 7 using phase portraits, Poincare maps 
and bifurcation diagram in the (/3 — X3) plane. The tran- 
sient hyperchaos exhibited by this circuit is describried 
in Sec. 8 and the hyper chaotic beats phenomenon is dis- 
cussed in Sec. 9. In Sec. 10 a multisim modelling of the 
memristive MLC circuit is made. With the help of the 
grapher facility of the software package, time series plots, 
phase portrait and power spectra which are qualitatively 
equivalent to the numerical observations are presented. 
In Sec. 11 a brief conclusion and possibilities of further 
studies are discussed. 



II. MEMRISTORS AND MEMRISTOR 
CHARACTERISTIC 

A Memristor is a two terminal circuit element supposed 
to have a functional relationship connecting charge q and 
fiux (j) and endowed with both a resistive property and a 
memory like nature. This memristor was first postulated 
by Leon Chua [28 as the missing fourth quartet in the 
fundamental circuit element family comprising of resis- 
tors, inductors and capacitors. Thirty seven years after 
Chua's prediction, researchers in the HP labs in 2008 re- 
vealed its existence by studying the electrical properties 
of certain nanoscale devices [29 . Though no combination 
of passive devices can reproduce the properties of any 
fundamental circuit element, much less that of a memris- 
tor, Chua nevertheless theorised that, in principle, mem- 
ristors with almost any desired (j) — q characteristic can 
be synthesised using analog circuits [28 . Based on this 
concept, Itoh and Chua [30] had proposed theoretically 
a three segment piecewise linear memristor with char- 
acteristic curves similar to those of Chua's diode, albeit 
with a difference. While the characteristic curve of the 
Chua's diode lies in the (v — i) plane, that of the proposed 
memristor falls in the {(t) — q) plane. Mathematically, this 
piecewise linear relationship is defined as 

q{cP) = Ga,cP + 0.5{Ga, - GaMl'P+M) - il'P-M)], (1) 

where Ga^ and are the slopes in the inner and outer 
segments of the characteristic curve. A typical character- 
istic curve of the Chua type memristor is shown in Fig. 
1(a). 

The current through this Chua-type fiux controlled 
memristor is given by 

i{t) = W{cP)v{t), (2) 

where W{(p)is the memductance of the memristor. It is 
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FIG. 1: (a) The characteristic curve of the memristor in the 
(0 ~ q) plane, (b)The variation of the memductance as 
a function of the flux (/)(t) across the memristor and (c) the 
variations of current and voltage across the memristor as a 
result of switching. 

so called because it has the units of conductance {U) and 
is given by 

where G^i , are as defined earlier. 

The memductance of the memristor takes on two par- 
ticular values as shown in Fig. 1(b), depending on the 
value of the flux across it. The higher memductance 
value can be referred to as the ON state and the 
lesser memductance value G^i as the OFF state. Ob- 
viously, as the flux across the memristor changes, the 
memristor switches or toggles between these two states, 
causing the current and voltage across the memristor to 
vary as shown in Fig. 1(c). This switching time vari- 
ant resistive property makes the memristor a desirable 
element in the modelling of nonlinear circuits. For ex- 
ample, memristor based chaotic circuits and their imple- 
mentation were described by [3TH34] . The coexistence of 
infinitely many stable periodic orbits and stable equilib- 
rium points has been reported in memristive oscillators 
by Messias and Botta [3F. Chaos and its control in a 
four dimensional memristor based circuit using a twin-T 
Notch filter has been studied in [36] . Transient chaos has 
been reported in a memristive canonaical Chua's circuit 
by [371 EH] • The dynamical behaviour and stability anal- 
ysis of a fractional-order memristor based Chua's circuit 
were made by Petras [39 . The nonlinear dynamics of 
a network of memristor oscillators were investigated by 
Corinto and Gilli [40 . Memristive chaotic circuits based 
on cellular nonlinear networks were studied by Buscarino 
et al [4T] . 

A. Multisim model of the Memristor 




FIG. 2: (a) Block diagram for a memristor based on basic 
principles and (b) a protypical analog circuit model of a three 
segment piecewise linear flux controlled active memristor used 
in the study. The parameters are Ri = lOKQ, R2 = lOOKQ, 
R3 = 100 KQ and C3 = 2.2nF for the integrator part, R4 = 
lOKQ for the output of the window comparator, R5 = 1450^1 
for the linear resistance and Rq = 1050^, R7 = 2KQ and 
i?8 = 2KQ for the negative conductance part. The reference 
voltages for the window comparator are fixed as ±1. 



resistance and a negative impedance converter (NIC) are 
switched ON and OFF alternatively based on the output 
pulse of a comparator. The comparator compares the 
flux through the memristor (that is the integral of the 
voltage across the memristor) between two reference lev- 
els, namely the breakdown points. For flux values lying 
within the upper and lower breakdown points (namely 
±1 flux units) the negative conductance is included in 
the circuit. For fiux values exceeding the breakdown 
points, the resultant of the linear resistance and nega- 
tive conductance, which are in parallel, is included in the 
circuit. The flux is obtained by an integrator circuit. By 
this action, the functional relationship between the flux 
and charge given in equation (1) is realized. A detailed 
description of this memristor model is given in [27] . 



While different models of the memristor haven 
been proposed by different authors, see for example, 
Muthuswamy and Kokate [32 , we make use of a protypi- 
cal analog circuit model of a three segment piecewise lin- 
ear flux controlled active memristor based on Faraday's 
law (j) — J^vdi designed by us, in an earlier work |27|. 
The block diagram of this memristor model, is given in 
Fig. 2(a). It is based on the time varying resistor (TVR) 
circuit proposed by Nishio and Mori [42 . Here a linear 



III. MEMRESTIVE 
MURALI-LAKSHMANAN-CHUA CIRCUIT 

A Memristive MLC circuit is constructed by replac- 
ing the Chua's diode by an active flux controlled mem- 
ristor as the nonlinear element in the classical Murali- 
Lakshmanan-Chua circuit. The schematic circuit is 
shown in Fig. 3(a), while the actual analog realization 
based on the prototype model for the memristor is shown 
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FIG. 3: (a) The memristive MLC circuit and (b) A Multisim 
Prototype Model of a memristive MLC circuit. The parame- 
ter values of the circuit are fixed as L = 21mH, R = 90011, 
C = 10. 5nF. The frequency of the external sinusoidal forc- 
ing is fixed 8.288kHz and the amplitude is fixed as 
F = 770mV {peak — to — peakvoltageVpp) . For the memristive 
part the parameters are as fixed earlier in Fig. 2. 



in Fig. 3(b). Applying Kirchoff's laws, the circuit equa- 
tions can be written as a set of autonomous ode's as 



dt 

'dt 
. di 

dp 
'dt 



C— = 



L— = 



i - W{(l))v, 

-V - Ri^ Fs'm{np), 

1. 



(4) 



where W{(j)) is the memductance of the memristor and is 
as defined earlier in Eq. ([3| and p is the time variable in 
the extended coordinate system. The normalized form of 
the circuit equations (7) are 



X2 = Xs - W{xi)x2, 

xs = -/3(x2 + /sin(a;x4), 
X4 — 1. 



(5) 



Here dot stands for differentiation with respect to the 
normalized time r (see below) and W{xi) is the normal- 
ized value of the memductance of the memristor, given 
as 



W{xi) 



dqixi) 
dxi 



ai, 

^2, 



Fi 
\xi 



< 1 
> 1 



(6) 



where a\ = Ga^/G and a2 = Ga^jG are the normalized 
values of Ga^ and Ga^ mentioned earlier and are negative. 
The rescaling parameters used for the normalization are 



C 



,X2 = V,X3 



,Xi 



Gp 



ui = 



G 
G~ 



,G 



1 



C 
2t:vC 

G 



C 



,(7) 



A. Stability Analysis 

When the driving force / = 0, the system can be con- 
sidered as a three dimensional autonomous system whose 
equilibrium state is given by 

A — {(xi, X2, ^3)1 X2 = Xs = 0, xi = constant}^ 

which corresponds to the xi-axis. The Jacobian matrix 
D at this equilibrium state is given by 



D, = 




(8) 



The characteristic equation associated with the system 
in this equilibrium state is 



cr^ -^P2(T +PlCr = 0, 



(9) 



where <j's are the eigenvalues and the coefficients p^'s are 
given as P2 = {P ca) and pi = /3(1 + a^). The eigen 
values that characterize the equilibrium states are given 
as 



CTl = U, (72,3 = ± 



4/3 



(10) 



where z = 1, 2. Depending on the eigen values, the nature 
of the equilibrium states differ. 

• When (/3— a^)^ = 4/3, the equilibrium state will be a 
stable/unstable star depending on whether (Z^ + a^) 
is positive or not. 

• When {P — > 4/3, the equilibrium state will be 
a saddle. muthulOb 

• When {P — < 4/3, the equilibrium state will be 
a stable/unstable focus. 



In this case the circuit admits self oscillations with nat- 
ral frequenc; 

v/(/3-a2)^-4/3 



1 r • • XI A/(/5-ai)^-4/3 

ural frequency varymg m the range — 2 < uJq < 



2 . It is at this range of frequency that the 
memristor switching also occurs. Therefore the memris- 
tor switching frequency uOmem depends on the normalized 
memductance values ai, a2. 



5 



IV. NON-SMOOTH SYSTEMS 

Dynamical systems that contain terms which are non- 
smooth functions of their arguments are cahed non- 
smooth or piecewise-smooth continuous systems. These 
systems are found to arise everywhere in nature, such 
as in electrical circuits that have switches, mechanical 
devices in which components (such as gear assemblies) 
impact with each other, mechanical systems with slid- 
ing, friction, etc. [^S", "IT. The common feature of these 
systems is that they are event driven, in the sense that 
the smoothness is lost at instantaneous events, for ex- 
ample upon the application of a switch. The phase 
space for such a piecewise-smooth system can be di- 
vided into two subspaces Si := {x G M^,x < 0} and 
5*2 := {x e W^^x > 0} by a discontinuity surface J^i2 
defined by I]i2 = x G : H{x) = 0. In each of the sub- 
spaces, the dynamics is governed by the smooth vector 
fields 



x{t) 



/i) if X G Si 
F2(x,/i) if X G 5'2. 



(11) 



where x denotes state variables. The surface 1112 is 
called the discontinuity set or discontinuity boundary or 
switching manifold. 

If 7^ F2(x,/i) at X = Si2, then we have the 

degree of smoothness to be one. Such systems with the 
degree of smoothness equal to one are called as Filippov 
systems [45j. However if = F2(x,/i) at x = T^i2 

but there is a difference in the Jacobian derivatives, such 
_^ dF2{x,^i) ^ _ Y,^^^ then the degree 



that 



dFi 



dx 



of smoothness is said to be two. Systems with de- 
gree of smoothness two or above are called piecewise- 
smooth flows. 

Generally as it is impossible to find explicit solutions 
to dynamical systems, it is convenient to introduce a flow 
function <I> generated by the vector field F{x^ ji) such that 
^{xo^t — to) corresponds to the point at time t on the tra- 
jectory that passes through Xq at time to- Let denote 
^ , where we have suppressed the parameter dependance 
for convenience. Then in terms of the flow function, a 
general dynamical sywhether grazing bifurcations induce 
a change in the dynamics from periodic to stem can be 
written as 



^t{x, t-to)= t - to)), ^(x, 0) 



(12) 



for all values of x and t. Then the unique solution for 
the dynamical system is given as x{t) = ^{xo^t — to). 

When a flow reaches a discontinuity surface, 
the vector field changes thereby causing a change 
in the flow function and a sudden qualitative 
change in the dynamics. These changes arise 
from a discontinuous jump in the Floquet multipli- 
ers and are called as Discontinuous Bifurcations or 
Discontinuity Induced Bifurcations (DIB)'s [44 . A 
grazing bifurcation is a particular form of discontinuity 



induced bifurcation which arises when a limit cycle of a 
flow becomes tangent to ( that is, it just grazes with) a 
discontinuity boundary. In the following sections we show 
that the memristive MLC circuit is a second order non- 
smooth system having two discontinuity boundaries and 
construct proper discontinuity mappings at these bound- 
aries to observe its dynamics. 



V. MEMRISTIVE MLC CIRCUIT AS A 
NON-SMOOTH SYSTEM 

The bistability nature of the memristor makes the 
memristive MLC circuit a piecewise-smooth or non- 
smooth continuous system of the second order. Referring 
to the memductance characteristic of the memristor, we 
find that the memristor switches states at x = +1 and 
at X = — 1 either from a more conductive ON state to 
a less conductive OFF state or vice versa depending on 
how we proceed across the discontinuity boundary. Con- 
sequently the phase space can be divided into three sub- 
spaces Si , 5*2 and Ss due to the presence of two switching 
manifolds Si ^2 and I]2,3 which are symmetric about the 
origin. These switching manifolds are defined by the zero 
sets of the smooth function i^i(x,/i) = (x — x*), = +1 
and H2 = (x — X*), X* = —1, respectively. The memris- 
tive MLC circuit can now be rewritten as a set of smooth 
ODE'S 



x{t) 



Fi(x,/i) if Hi{x,ii) > and H2{x,ii) < 



F2{x,ii) if Hi{x,ii) < and H2{x,/i) > . 

(13) 

where ja denotes parameter dependance of the vector 
fields and the scalar functions. The vector fields F^'s 
are 



F^(x,/i) = 



X2 

xs - aiX2 
{x2 + Xs) + fsin{(jJX4) 
1 



,i=l,2. (14) 



The conditions for grazing of the periodic orbits (re- 
fer Eq. (26) in the Appendix A) are satisfied at both 
the discontinuity boundaries 111,2 and S2,3- Hence one 
can check that grazing bifurcations indeed occur in the 
memristive MLC circuit at both of these discontinuities. 

Further, the switching action of the memristor induces 
a tendency in the system to exhibit self-oscillations. Such 
self-oscillations are known to occur in nonsmooth me- 
chanical systems and were studied extensively by [46] . 
Even the frequencies and amplitudes of such oscillations 
were estimated using specialized methods. In this sys- 
tem described by Eq. ([5|, if an external periodic forcing 
is applied, there may occur an interaction between this 
externally induced oscillations and the self-oscillations re- 
sulting in the occurrence of beats in the system. This has 
actually been observed and is described in detail in Sec. 
9. 
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A. Discontinuity Mappings 



where the PDM correction 7 is given as 



Generahy non-smooth discontinuity induced birfuca- 
tions are difficult to analyse in piecewise-smooth systems, 
because it is necessary to establish the fate of topolog- 
ically distinct trajectories close to the structurally un- 
stable event that determines the bifurcation. To over- 
come this difficulty, the concept of the discontinuity map 
was introduced by Nordmark [47]. This is a synthe- 
sised Poincare Map whose orbits completely describe the 
dynamics of the system. Depending on the change in 
time or otherwise, one can construct two discontinuity 
maps, namely the Zero-time Discontinuity Map (ZDM) 
and Poincare Discontinuity Map (PDM). For details of 
deriving anallytically these mappings, one may refer to 
di Bernado et al. [44 . 

If we consider the grazing point as Xq = 
{^10,^:20,^30,^40} = {±l,0,X3o,X4o} then the vector 
fields of the system become continuous at the grazing 
point, that is, Fi{x) = F2{x) = F{x) where 



/ 



Fix) 



I 



\ 




(15) 



Therefore the system is continuous at the grazing point, 
but is discontinuous in the first derivative, that is 

dx dx 

This makes the memristive MLC circuit a second order 
non-smooth system. Following di Bernado et al. [44 , the 
ZDM and PDM for the memristive MLC circuit can be 
derived (for details see Appendix A). The Zero Disconti- 
nuity Map is derived as 



X ^ ZDM 




Here 



Hrnin (^) 



dH 

dx 



if Hmin{x) > 
if Hmin{x) < . 



x^O{\xf) 



(16) 



where 0(|xp is the correction term in the Taylor expan- 
sion of H{x). The ZDM correction 5 is given as 



ai 



■ a2)(ij)'/'(ai - 2(a2 + P))x3 + 2fsin{ujXi 



\ 







(17) 

Obviously we find that the correction to the time variable 
X4 is zero, thereby justifying the nomeclature ZDM. In a 
similar manner the Poincare Discontinuity Map is derived 
as 



^ PDM = 



X + 7 



if H{x) > 
if H{x) < . 



(18) 



2(ai 



\ 

-a2)(Qi+/3) ^/xT 



\ 








(19) 



It is essential that these corrections are to be incorpo- 
rated while deriving the analytical solutions or while per- 
forming the numerical integration, so as to obtain a true 
picture of the non-smooth dynamics. 



B. Event Driven Numerical Simulations 

As the system under consideration is a low dimensional 
system with just two discontinuity boundaries, disconti- 
nuity corrections similar to the above may also be made 
by following direct numerical simulation using what are 
known as event driven schemes [44 . Often these meth- 
ods are fast, realistic and accurate. Here the trajec- 
tories in the regions Si are solved using the standard 
numerical integration algorithms for smooth dynamical 
systems such as the Runge-Kutta ODE solver algorithm. 
The discontinuity boundary defined as a zero set of a 
smooth function Hij{x) = is identified and the inte- 
gration is switched upon crossing the boundary to follow 
the smooth dynamics pertaining to that particular sub- 
space and appropriate corrections, which are determined 
by forward and backward integrations, are applied. This 
scheme therefore reduces the time-integration of a tra- 
jectory of a piecewise smooth dynamical system to the 
problem of finding a set of event times tk and events 
h\^^ such that 

To implement this one sets up a series of monitor func- 
tions^ the values of which are computed during each in- 
tegration step. When one of these monitor functions 
change sign, a root finding method is used to find ac- 
curately where Hij = 0. The corrections are then deter- 
mined by logical manipulations. In our investigations, we 
have used both the above methods to validate the concept 
of discontinuity corrections as applied to the memristive 
MLC circuit. 



VI. GRAZING BIFURCATIONS 

In the presence of forcing, we find that the circuit ad- 
mits limit cycle motion resulting in periodic orbits. When 
a periodic orbit impacts the discontinuity surfaces tan- 
gentially, a local stretching of the phase space leading to 
a profound destabilizing effect on the dynamics have been 
observed. We find that the periodic orbit after impact 
is thrown away from the discontinuity surfaces to regions 
infinitely deep in the subspaces or 5*3. Hence recording 
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FIG. 4: Grazing bifurcations for the system defined by Eq. 



14) with ZDM and PDM corrections given by Eqs. & 
19} for the parameters ai = -0.465116, a2 = -0.120289, 
/3 = 0.163613, uj = 0.325940 and / = 0.747114. A peri- 
odic orbit (a) for the initial condition {7.0, 7.0, 7.0, 0.0} as 
it approaches the grazing boundary S = 1 flux units in the 
subspace S2,. (c) at the onset of grazing for the initial con- 
dition {4.5,4.5,4.5,0.0} and (e) after grazing for the initial 
condition {1.5,2.5,2.5,0.0}. The Lyapunov exponents A-s 
in each case are {0.000000, -0.022011, -0.022012,0.000000}, 
{0.000040,-0.0218890,-0.021894,0.000000} and {0.000281, 
-0.021557, -0.021561,0.000000} respectively. Figs. (b),(d) 
and (f) are the blown up portions of the orbit in Figs, (a), 
(c) and (e) respectively. The insets in Figs, (d) & (f) show 
the trajectories at grazing and after complete onset of grazing 
with more clarity. 



this grazing bifurcation numerically is difficult. To over- 
come this problem, the parameters are fixed at some de- 
sired values and the initial conditions are changed start- 
ing from faraway regions in subspaces or S^, , to regions 
close to the discontinuity boundaries. By following this 
approach it is possible to find the behaviour of the peri- 
odic orbit in phase space, before impact, just at the in- 
stant of impact and after the impact has taken place with 
the discontinuity boundary. The parameter values are 
fixed as ai = -0.465116, as = -0.120289, f3 = 0.163613, 
UJ = 0.325940 and / = 0.747114. For an initial condition 
Xq — {7.0,7.0,7.0,0.0} we find a periodic orbit close to 
the discontinuity boundary. This is shown in Fig. 4(a). 
A blown up portion of a section of the periodic orbit 
is shown in Fig. 4(b) for clarity. The Lyapunov ex- 
ponents for this periodic orbit are Ai = 0.000000, A2 = 
-0.022022, A3 = -0.022012 and A4 = 0.000000. As the 
initial point is brought closer to the discontinuity bound- 
ary, grazing sets in. This is shown in Fig. 4(c) & 4(d) for 



the initial condition Xq = {4.25, 4.25, 4.25, 0.0}. The cor- 
responding Lyapunov exponents are Ai = 0.000040, A2 = 
-0.021890, A3 = -0.021894 and A4 = 0.000000. For an 
initial condition Xq = {1.5, 2.5, 2.5, 0.0} we find that the 
periodic orbit has been blown up into a fully chaotic 
attractor with Ai = 0.000281, A2 = -0.021557, A3 = 
-0.021561 and A4 = 0.000000. These are shown in Figs. 
4(e) & 4(f). Thus we find that the attractor which is 
purely periodic before transversal intersection with the 
discontinuity boundary, becomes chaotic as a result of 
grazing bifurcations. This result is also proved by apply- 
ing the K-test for chaos as described in the next section. 



A. 0-1 Test for Chaos 

Recently Gottwald and Melbourne [48 have intro- 
duced a new test called the binary test or the 0-1 test 
or the K-test for chaos. This test uses a scalar time se- 
ries data obtained from the system to derive a rotation 
coordinate d and translation coordinates (p, q) to form a 
Euclidean group E(2\ As reported by Nicol et al. 
the group extensions of this Euclidean group are sublin- 
ear and are 

• bounded when the system dynamics is periodic or 
quasiperiodic and 

• unbounded random Brownian motion like be- 
haviour when the system dynamics is chaotic 

This causes the mean square displacement K of points 
in the space of the translation variables, namely the 
{p — q) space, to become when the system behaviour is 
periodic and 1 if it is chaotic. The main advantages of 
this test are that 

• the dimension of the dynamical system and the 
form of the underlying equations are not important 

• the determination of embedding dimension and de- 
lay time, the essential requisites of phase space re- 
construction, are done away with 

• there is an element of definiteness. This is because 
a value of 0.01 for K denotes definitely a nonchaoic 
behaviour, which is not possible with the lyapunov 
exponents 

• it provides for a straight forward visualisation of 
the dynamics in the translation variables {p — q) 
space. 

The success of this procedure depends on leaving long 
transients such that the system behaviour settles down 
to its asymptotic state and having larger values of in- 
tegration time N. The test is universally applicable to 
any deterministic dynamical system, in particular to or- 
dinary and partial differential equations, and to maps 
[48], as well as to deterministic systems with noise j5Q] . 
This test has also been applied successfully to Chua's 
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circuit [ST, to fractional Chua's citcuit [52j [53] and to 
fractional simplified Lorenz oscillator [54]. 

Due to the advantages mentioned above, the 0-1 test 
is getting greater acceptance by nonlinear researchers for 
characterizing the dynamics of nonlinear systems. 



grazing and K = 0.000661 for behaviour at grazing and 
K = 0.998337 for the system's behaviour after grazing. 
These values of K prove that grazing bifurcations occur- 
ring in the system are indeed responsible for the change 
in the dynamics of the system from a periodic behaviour 
to chaotic behaviour. 



B. Confirmation of Grazing Bifurcation Induced 
Dynamics using 0-1 Test for Choas 

The 0-1 test can be applied to the memristive MLC 
circuit described by Eq. (14) along with the ZDM and 
PDM corrections given by Eqs. (17) & (19), to find 
out whether grazing bifurcations induce a change in the 
dynamics from periodic to chaotic behaviour, indepen- 
dent of the values of the associated Lyapunov exponents. 
To apply this test, the translation coordinates p{t)hq{t) 
given by [48 are obtained as 



Pc(n) = ^0(j)cos(i9c), 



(lc(n) = ^^(j)sin(l9c), 



OcU) = jc + ^(/>(z). 



(20) 



where 0(j) is the scalar time series for any variable of the 

system of length j = 1,2,3 c G (0,7r). The mean 

square displacement of the translation variables is given 
for a given value of c as 



Me 



^ N-n 

c(n) = lim — — - V [Pc(j+n)+pc(j)]^+fec(j+n)+^c(j)]^ 



(21) 

where n = 0,1,2 N/10. The above mean square 

displacement grows linearly for chaotic behaviour but re- 
mains bounded for nonchaotic behaviour. The asymp- 
totic growth rate of this quantity can be found from lin- 
ear regression as 



Kr = lim 



logMcjn) 
log{n) 



(22) 



The median of Kc for c G (0, tt) denoted as K gives the 
desired characterization parameter. 

The attractors of the memristive MLC circuit in the 
space of the translation variables, namely the {p — q) 
space for the periodic behaviour before grazing, for the 
behaviour at grazing and for the chaotic behaviour af- 
ter grazing are shown in Figs. 5(a), 5(c) & 5(e) respec- 
tively. The corresponding variations in the values of the 
mean square displacement Kc as a function of the con- 
stant c are shown in Figs. 5(b), 5(d) & 5(f). The median 
K of these Kc values is obtained for each case and are 
given 8iS K = 0.000440 for the system's behaviour before 
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FIG. 5: Attractors of the system defined by Eq. (14) in the 



plane of the translation variables, namely the (p — q) plane 
for (a) periodic behaviour before grazing , (c) for the peri- 
odic behaviour at grazing and (e) for the chaotic behaviour 
after grazing. (b),(d) & (f) the corresponding variations in 
the mean square displacement K values as a function of the 
constant c. 



VII. CODIMENSION-M BIFURCATION 

A codimension-m bifurcation refers to a bifurcation 
scenario induced by the concurrent variation of m pa- 
rameters of the system. To unfold this bifurcation fully 
one requies dim — dimensional parameter space [55] . Let 
us consider the case where the the values of the circuit el- 
ements are fixed as C = 12.5nF, L = 19.1mi^, = 27rz^, 
where u = 8300Hz is the frequency and F = 1000 mVpp 
the peak-to-peak value of the external force. As the series 
resistance R enters into the expressions of all the normal- 
ized parameters given in Eq. (7), it is but natural that 
all of them change when the resistance alone is varied. 

For example, when R = 720^1, the normalized param- 
eters of the system are ai = —0.669767, a2 = —0.173216, 
/3 = 0.339267, u = 0.469354, / = 0.119949 and h = 
0.013387. For this choice of parameters, the behaviour 
is hyper chaotic with the associated Lyapunov expo- 
nent values as (Ai = 0.136691, A2 = 0.034005, A3 = 
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FIG. 6: (a) The phase portrait of the hyperchaotic attractor 
in the (x2 — X3) plane and (b) its corresponding Poincare 
map. The series resistance for this case is R — 720Q giving 
the normahzed parameter values as ai = —0.669767, a2 = 
-0.173216, 13 = 0.339267, uj = 0.469354, / = 0.119949 and 
h = 0.013387. 



-0.109690, A4 = 0.000000). The phase portrait of this 
hyperchaotic attractor in the {x2 —^3) plane and its cor- 
responding Poincare map are shown in Figs. 6(a) & (b). 




FIG. 8: The power spectrum of (a) the hyperchaotic attractor 
and (b) the periodic attractor. 
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FIG. 7: (a) The phase portrait of the periodic attractor in 
the (x2 — X3) plane and (b) its corresponding Poincare map. 
The series resistance for this case is = 788Q causing the nor- 
malized parameters to be ai — —0.733002, a2 = —0.189577, 
13 = 0.406377, uj = 0.513682, / = 0.143676 and h = 0.012232. 

However when R = 7881^, the normalized parameters 
get changed as ai = -0.733002, as = -0.189577, /3 = 
0.406377, u = 0.513682, / = 0.143676 and h = 0.012232. 
The phase portrait of the system for this choice of param- 
eters in the {x2—xs) plane and its corresponding Poincare 
map are shown in Figs. 7(a) & (b). The values of the as- 
sociated Lyapunov exponents are (Ai = —0.000052, A2 = 
-0.116821, A3 = -0.117559, A4 = 0.000000). 

The power spectrum for the system shown in Fig. 8(a) 
when the series resistance R = 720Q is of broadband na- 
ture, indicating the hyperchaotic nature of the attractor. 
However the power spectrum in Fig. 8(b) for the case 
when R = 788^1, contains a single peak, which clearly 
reveals the behaviour of the system to be periodic. 

Further the hyperchaotic and the periodic natures 
of the system's behaviour are also confirmed by the 
phase portraits in the plane of the translational variables, 
namely the {p — q) plane as shown by Figs. 9(a) & (c) 
and in the variation of the K values shown in Figs. 9(b) 
& (d), for a constant value of the random constant c in 
Eq(20), for a range of values of the discrete time N. It 



has been found that the median of these K values for 
the hyperchaotic attractor is = 1.04475 and for the 
periodic case K = 0.003138. 



FIG. 9: The phase portraits of (a) the periodic attractor and 
(c) the hyperchaotic attractor in the plane of the translational 
variables, namely the {p — q) plane, (b) & (d) The variation 
of the K values as a function of discrete time N for these two 
cases respectively. 



The stable hyperchaos and periodicity exhibited by the 
system is a consequence of the nonsmooth grazing bifur- 
cation undergone by the system. 



A. Bifurcation Diagram 

As it is difficult to depict this codimension-5 bifur- 
cation in three dimensional space, we capture this bi- 
furcation with the help of a bifurcation diagram in a 
two dimensional plane. For this purpose we vary the 
series resistance in the system in Fig. 3 in the range 
400(7 < R < 100017, starting from 100017. Due to this 
variation in resistance, the normalized parameters vary 
as follows, -0.930233 < ai < -0.372093, -0.240577 < 
a2 < -0.096231, 0.104712 < (3 < 0.654450, 0.260752 < 
Uext < 0.651880 and 0.037021 < / < 0.231383. These 
variations in parameters cause the system to transit from 
a stable periodic attractor to a hyper chaotic attractor 
suddenly at a resistance value R = 74017 or equivently 
/3 = 0.358377 . Hence we can say that the system un- 
dergoes a codimension-5 bifurcation. The circuit then 
alternates between hyper chaotic, chaotic and transient 
hyper chaotic states as shown in the in the {P — xi) plane 
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representation of the said bifurcation in Fig. 10(a). The 
transient hyper chaotic states are shown by blank re- 
gions or gaps in the bifurcation diagram. This is because 
at the instants of transition from hyper chaotic to pe- 
riodic state, the attractor is thrown from the regions in 
subspaces Si, S2 and close to the fixed point at the 
origin to infinitely faraway regions in either the subspace 
or where it eventually settles into a chaotic state. 
The Lyapunov spectrum for the same range as in the bi- 
furcation diagram is shown in Fig. 10(b). This spectrum 
follows closely the dynamics represented by the bifurca- 
tion diagram. 
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the dynamics to remain the same. Also the results ob- 
tained using these step sizes are similar to those one may 
get using the standard value for the step size, namely 
h = 0.01. 





5.0 



FIG. 11: (a) The time series plot showing transient hyper- 
chaos and (b) the Local Lyapunov spectrum confirming the 
same. The parameters cosidered are ai = —0.725581, a2 = 
-0.187650, 13 = 0.398168, ujext = 0.505404, / = 0.200602 and 
h = 0.012432. 



FIG. 10: (a) The bifurcation diagram in the {/3 — X3) plane 
corresponding to the codimension-5 bifurcation and (b) the 
corresponding Lyapunov spectrum. 

In order to ensure that our numerical simulation fol- 
lows closely the dynamics, we have taken the step size as 
h = {27t)/{ujN), where we have chosen N = 1000. This 
gives the value of the step size as h = 0.0193. Also as the 
system represented by Eq. (5) is a non-autonomous sys- 
tem, this choice of step size is not altogether new, as this 
choice has to be necessarily made if one wishes to take 
a stroboscopic section of the attractor, see for example 
|56J. However the dynamics is independent of the step 
size of the numerical integration. To verify this, we have 
also repeated our simulations with a smaller step size, 
by choosing N = 1500 and N = 2000 and have found 



VIII. TRANSIENT HYPERCHAOS 

Transient chaos is a ubiquitous phenomenon observed 
in many nonlinear dynamical systems, wherein a trajec- 
tory behaves chaotically for a finite amount of time before 
settling into a final (usually nonchaotic) state. This arises 
due to the presence of nonattracting chaotic saddles in 
phase space [57-62 . It is known that chaotic saddles 
and transient chaos are responsible for important physi- 
cal phenomena such as chaotic scattering [63j, and par- 
ticle transport in open hydrodynamical flows [64 . They 
are belived to be the culprit for catastrophic phenomena 
such as voltage collapse in electrical power systems [65] - 
67 and species extinction in ecology [68]. An extensive 
study of chaotic transients in spatially extended systems 
was studied by Tel and Lai [69^ . It has also been reported 
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FIG. 12: The phase portrait in the (xi — X2) plane for (a) 
the transient hyperchaotic state upto approximately 2.2 x 10^ 
normalized time units and (b) in the asymptotic state beyond 
2.2 X 10^ normalized time units. 



solution with greater accuracies, the attractor continues 
to fill out, thus justifying the low positive values of the 
largest local Lyapunov exponent in the asymptotic state. 
This is a consequence of the non-smooth nature of the 
system. 
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in a memristive Chua oscillator by [371 EH] • In our present 
study we have found transient hyperchaos in the memris- 
tive MLC circuit. To observe this we set the parameters 
of the circuit as ai = —0.725581, a2 = —0.187650 and 
/3 = 0.398168. The frequency of the external force is 
set as ujext = 0.505404 and the strength of the force as 
/ = 0.200602. The step size for the numerical integration 
is assumed as = 0.012432. The initial conditions are 
chosen as Xq = {1.0,-0.1,0.1,0}. For these values of 
the parameters, the circuit initially exhibits hyperchaos 
and then switches over to chaotic oscillations. This can 
be seen in the time series plot shown in Fig. 11(a). The 
local Lyapunov spectrum, drawn as a function of time 
shown in Fig. 11(b), illustrates clearly this transition. 
While the two largest local Lyapunov exponents (LLE) 
are positive initially, at the instant of the transition to 
chaos, the second largest local exponent changes from 
positive to negative values. The phase portraits of the 
oscillator in the Xi — X2 plane in the transient hyper 
chaotic and asymptotically chaotic states are shown in 
Figs. 12(a) & 12(b). Though the attractor shown in Fig. 
12(b) appears as a period 1-T orbit, it is in fact an over- 
lay of different trajectories. As one tries to compute the 



FIG. 13: (a) The time series plot showing persistent or sta- 
ble hyperchaos for a slightly different set of parameters ai = 
-0.725581, as = -0.187650, /3 = 0.398168, ojext = 0.508467, 
/ = 0.183006 and h = 0.012357 and (b) the Local Lyapunov 
spectrum confirming the same. 

In contrast to the transient hyperchaos, persistent 
or stable hyperchaos observed for slightly different pa- 
rameter values ai = —0.725581, a2 = —0.187650 and 
13 = 0.398168, ujext = 0.508467,/ = 0.183006 and 
h — 0.012432 is shown by the corresponding time se- 
ries, local Lyapunov spectrum in Figs. 13(a) & 13(b) 
and phase portraits in Figs. 14(a) & 14(b). 



IX. HYPERCHAOTIC BEATS 

Beats is a phenomenon arising out of the interference 
between two different oscillations having small differ- 
ence in frequencies. This phenomenon has been widely 
studied in linear systems. However chaotic and hyper- 
chaotic beats are found to occur in nonlinear systems 
also. Chaotic and hyper-chaotic beats were identified for 
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FIG. 14: The phase portrait in the (xi — X2) plane for (a) the 
persistent hyperchaotic state at approximately 2.2 x 10^ nor- 
malized time units and (b) in the asymptotic region beyond 
2.2 X 10^ normalized time units. 



sitor are shown in Figs. 15(a-b). The power spectra of 
xi variable as well as that of the memristor switching 
are shown in Figs. 15(c-d), respectively. The time series 
of the modulated normalized flux xi is shown for an ex- 
tended range in Fig. 16(a), while a section of the same is 
shown in Fig. 16(b). Similar results obtained by numer- 
ical simulation employing even- driven scheme are shown 
in Figs. (17 & 18). 
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FIG. 15: Hyperchaotic Beats: Numerical Simulations with 
Analytical Corrections : The phase portrait in (a) the 
(xi — X2) plane, (b) the variation in time of the mem- 
ductance VK(xi), (c) the power spectrum of the normal- 
ized flux xi across the memristor and (d) the power spec- 
trum of memristor switching showing it to be chaotically 
time varying resistor (CTVR). The parameters are /3 = 
0.355263, ai = -0.697674, as = -0.180433, ujext = 0.508938 
and / = 0.081643. The Lyapunov exponents A^s are 
{0.125510, 0.016904, -0.017554, 0.000000}. 



the first time in a system of coupled Kerr oscillators and 
coupled Duffing oscillators with small nonlinearities and 
strong external pumping [70 . Chaotic beats were also re- 
ported in coupled non-autonomous Chua's circuits [71 . 
Weakly chaotic and hyper-chaotic beats were reported in 
individual and coupled nonlinear optical subsystems, re- 
spectively, describing second harmonic generation (SHG) 
of light by Sliwa et al [72] . In all these coupled systems 
[7QH72] the occurrence of beats is explained as due to the 
interaction between the self oscillations or driven oscilla- 
tions of each of the coupled subsystems. Using extensive 
Pspice simulations, the occurrence of beats in individ- 
ual driven Chua's circuit with two external excitations 
has been reported [TSj |74] . Chaotic beats have also been 
found to occur in individual Murali- Lakshmanan-Chua 
(MLC) circuit by the same authors [17]. In these cases of 
single systems [TTl [73l [74] , two different external driving 
forces with slightly differing frequencies were the cause 
for the occurrence of beats. Chaotic beats in a Memris- 
tive Chua's circuit also was reported earlier by [27 . In 
the present work, we have identified hyper chaotic beats 
in the memristive MLC circuit. The phase portrait in 
the {xi — X2) plane, the chaotic switching of the memr- 



The presence of beats in any electronic circuit is ex- 
plained in detail by Nikolai [75 as due to the incomplete 
synchronisation between self oscillations and the exter- 
nal forcing caused by a weak perturbation. The mem- 
ristive MLC circuit can be considered to contain a res- 
onant part, comprising of the L — C combination and 
the memristor. In addition to this it has a nonlinear 
part whose role is solely taken by the memristor itself. 
The frequency of the circuit is determined by the reso- 
nant combination. The dissipative nature of the resonant 
part tries to destablise the stationary state of the circuit, 
while the non linearity, which is basically a time vary- 
ing negative conductance is an active element and hence 
tries to restore the amplitude. Due to this interaction, a 
particular mode of periodic oscillations having a band of 
frequencies is preferred. Due to the existence of a band 
of resonant frequencies, the characteristic frequencies of 
this mode will change with a small variation of parame- 
ters of the circuit. An external periodic force applied to 
the circuit will modulate the environment in which the 
intrinsic mode of oscillations lives. As a result of this 
modulation, the non linearity and dissipation will again 
tend to change the intrinsic mode of oscillation. If the 
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FIG. 16: Hyperchaotic Beats: Numerical Simulations with 
Analytical Corrections : (a) The time series of the hyperchaot- 
ically modulated X2 variable and (b) an expanded portion of 
the same for a short stretch of time. 



perturbations introduced by the external periodic force 
is insufficient for proper change of intrinsic mode, then 
the circuit will be in and out of synchronisation inter- 
mittently and will demonstrate a partial collapse and re- 
vival of amplitudes of the circuit variables which we call 
as beats. While the presence of beats may be attributed 
due to this interaction between the self oscillations of the 
circuit and the external forcing, the nonsmooth nature of 
the circuit and the presence of grazing bifurcation of peri- 
odic orbits at the topological discontinuities are believed 
to endow the beats phenomenon with a hyper chaotic 
nature. Nevertheless, purely quasiperiodic beats may be 
observed when the initial conditions are so chosen that no 
grazing bifurcations occur. To observe the hyper chaotic 
beats numerically we set the parameters as /3 = 0.355263, 
ai = -0.697674, as = -0.180433. The eigen values of 
the unforced case (Sec. 3.1) for this choice of parameters 
are 

cTi = 0.0, |xi I < 1, (72,3 = -0.526469 ± i0.570922, 

0-1 = 0.0, \X2\ > 1, 0-2,3 = -0.267848 ± zO.589594. (23) 




FIG. 17: Hyperchaotic Beats: Numerical Simulations (Event 
Driven Scheme) : The phase portrait in (a) the (xi — X2) 
plane, (b) the variation in time of the memductance VK(xi), 
(c) the power spectrum of the normalized flux xi across the 
memristor and (d) the power spectrum of memristor switching 
showing it to be chaotically time varying resistor (CTVR). 
The parameters are p — 0.355263, ai = —0.697674, a2 = 
-0.180433, i^ext = 0.528730, / = 0.081643 and h = 0.011884. 

Here both the sets 

M = {{xi,X2,X3)\ \xi\ < 1,X2 = X3 = 0}, 

and 

^1,3 = {(^1,^2,^3)1 1^1 1 > l,X2 = X3 = 0}. 

are spirals or the system exhibits damped self oscillations 
with a natural frequency uJq which is found to vary in the 
range 0.156187 < Uq < 0.567041. It is at this frequency 
that the memristor switching also occurs. This means 
^mem = ^o- When an external driving force is applied, 
we find that the unstable damped oscillations stabilize 
to limit cycles and when the frequency of the external 
source is adjusted close to the memristor frequency, the 
phenomenon of beats occurs. Here in this circuit, when 
the normalized frequency of the external source is fixed 
as Uext = 0.508939, and the normalized amplitude of 
the driving force is kept as / = 0.057846, while all the 
other parameters are kept fixed as above, hyper chaotic 
beats occur with the Lyapunov exponent values (Ai = 
+0.125510, A2 = +0.016904, A3 = -0.107554 and A4 = 
0.000000) as shown in Figs. (15) & (16). 

X. EXPERIMENTAL STUDY 

We present here the results of the circuit simulation of 
hyper chaotic beats in the memristive MLC circuit us- 
ing the Multisim model of the memristor mentioned ear- 
lier. The circuit implementation itself is shown in Fig. 
3(b). The parameter values of the circuit are fixed as 
L 21mH, R = 900n, C 10.5nF. For the mem- 
ristor, the parameters are Ri = lOKQ^ R2 = lOOKQ, 
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FIG. 18: Hyperchaotic Beats: Numerical Simulations (Event 
Driven Scheme) : (a) The time series of the hyperchaotically 
modulated X2 variable and (b) an expanded portion of the 
same for a short stretch of time. 



Rs = 100 KQ and Cs = 2.2nF for the integrator part, 
i?4 = lOKQ. for the output of the window comparator, 
= 1450^1 for the linear resistance and Rq = 1050^1, 
Rj = 2KQ and Rs = 2KQ for the negative conductance 
part. The reference voltages for the window compara- 
tor are fixed as ilF. The frequency of the external si- 
nusoidal forcing is fixed as Uext = 8.288KHz and the 
amplitude is fixed as F = 770mVVpp. For these choice 
of parameters, we find hyper chaotic beats. The phase 
portrait in the {(j) — Vc) plane and the time series plot 
of the memristor switching are shown in Figs. 19(a-b). 
The power spectrum of the flux (j) across the memristor 
as well as that of the memristor switching are shown in 
Fig. 19(c-d) respectively. The time series of the volt- 
age across the capacitor for an extended range is shown 
in Fig. 20(a) and for a short span in Fig. 20(b). All 
these are similar to the numerical simulation results of 
the hyper chaotic beats phenomenon reported in Sec. 9. 



XI. CONCLUSION 

In this study some of the questions that were consid- 
ered at the begining have been addressed. Firstly the 
prototypical multisim model of the memristor, which we 
had proposed earlier and even designed a nonautonomous 




FIG. 19: Hyperchaotic Beats: Multisim Simulations : The 
phase portrait in (a) the {(j) — vc) plane, (b) the variation in 
time of the memductance W {(/)), (c) the power spectrum of 
the flux (j) across the memristor and (d) the power spectrum of 
memristor switching showing it to be chaotically time varying 
resistor (CTVR). 



Chua's memrestive oscillator using it, has been success- 
fully employed in designing a memristive MLC circuit. 
We have identified this memristive MLC circuit as a 
piecewise smooth system of second order with two dis- 
continuous boundaries. Grazing bifurcations which are 
a form of non-smooth bifurcations, have been found to 
occur in this circuit giving rise to chaotic, hyper chaotic 
and transient hyper chaotic dynamics. The bist ability 
nature of the memristor, the interactions of the self os- 
cillations with the external forcing and grazing bifurca- 
tions have given rise to hyper chaotic beats. These have 
been verified by numerical simulations after constructing 
proper ZDM and PDM mappings. Multisim modelling 
of the circuit has resulted in the verification of the con- 
cepts proposed and behaviors observed through numeri- 
cal computations. 

This work is important in that it has brought the in- 
vestigation of analog circuits within the framework of 
non-smooth bifurcation theory. Earlier such bifurcations 
were identified and studied in detail in mechanical sys- 
tems or power electronic circuits alone. As memristor 
has an inherent bist ability, all memristive systems may 
supposed to be endowed with a non-smooth nature, and 
an analysis of these may help in unravelling their true 
dynamics, which would not be possible when one applies 
the theories and concepts developed for smooth systems. 
Already we have performed this analysis for the mem- 
ristive Chua's oscillator and driven memristive Chua's 
oscillator. The results of this study will be published 
later. 

As further study, it is clear that the insights gained by 
the analysis of individual and networks of such oscillators 
may help in the building up of potential applications such 
as neuromorphic circuits and dense nonvolatile and dy- 
namic memories. Further we feel that the identification 
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FIG. 20: Hyperchaotic Beats: Multisim Simulations : (a) The 
time series of the hyperchaotically modulated vc variable and 
(b) an expanded portion of the same for a short stretch of 
time. 



of the hyper chaotic saddles and the control of transient 
phenomenon as well as use of hyper chaotic beats to phys- 
ical and real time applications are some of the areas that 
can be further investigated. 
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Appendix: ZDM and PDM Corrections for Grazing 
Bifurcations 



In this section the conditions for grazing bifurcations 
to occur in a general second order non-smooth system 
having a single discontinuity boundary and the analyt- 
ical expressions for the discontinuity mappings [44] for 
transversal crossings of the periodic orbits are given. Us- 
ing these, it is shown that, the ZDM correction and PDM 
correction given in Eqs. (16) & (17) respectively for the 



Let us consider a non-smooth system of n*^ order, hav- 
ing a single discontinuity surface defined by 



S := {x G D : H{x) = 0}. 



(24) 



In the neighbourhood of this boundary the non-smooth 
system can be written locally as a system of smooth 
ODE'S 



\f2{x,) iiH,{x)<0. 



(25) 



For the dynamical system under consideration, the vec- 
tor field Fi{x) generates a flow function while the vec- 
tor field F2 will generate a flow ^2- Hence to detect this 
cross-over, we assume a scalar smooth function H{x)^ 
such that the zero set of this smooth scalar function de- 
fines the discontinuity surface E12. A point x = is said 
to be a grazing point of a flow if the three conditions that 
follow are satisfied. 



v{x*) 
a{x*) 



0, 

d_ 
dt 

dH{x 



H{<^>iix*,0) 

Fi{x*), 



dx 

dH{x) dFijx*) 
dx dx 
:= a* > 0. 



(26) 



All of these expressions are evaluated at the point x = x^ . 

The first condition states that G S, whereas the 
second condition states that the flow is tangent to S at 
X*. The third condition gives the sense of the direction 
of the flows upon reaching the discontinuity surface and 
crossing it. A positive value, a(x*) > 0, tells that the flow 
is convex to the boundary. However, a negative value for 
a(x*) < 0, tells that the flow is concave to the boundary. 
In addition to these there is a transversality condition 
that is to be satisfied, namely 



X 



dx 



Fi{x) 



OH^jx) 
dx 



F2{x) > 0. 



(27) 



memristive MLC circuit can be obtained. 



This condition ensures that no sliding motion occurs 
along the discontinuity boundary. 

The unique solution of the dynamical system can be 
obtained by monitoring the scalar function H{x) during 
the numerical integration of the system equations and 
when H{x) = 0, the so called Zero Discontinuity Map 
(ZDM) is applied. The zero discontinuity mapping takes 
into account both the changes in the state variables and 
the vector field before and after the discontinuity jump. 

To obtain a periodic orbit we can reason out that it is 
possible to find at least one point, say which lies not 
on the discontinuity boundary but on a local transversal 
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Poincare surface such that x* G Hp. Here the lo- 
cal transversal Poincare surface is defined using a scalar 
function Hj^pix^) such that 



(28) 



Solving these equations using Newton's method, the pe- 
riodic orbits can be located. The periodic motions with 
grazing impacts on the discontinuity surface will then be- 
come a fixed point of this local transversal Poincare map- 
ping. These motions are sensitive to even small pertur- 
bations leading to eventual bifurcations. By constructing 
proper discontinuity maps, these grazing bifurcations can 
be simulated numerically. 

Let X* be a regular grazing point of a second order 
piecewise-smooth system defined by Eq. (25). Then the 



Zero Discontinuity Map describing the trajectories of this 
system in a neighbourhood of the grazing trajectory has 
a 3/2-type singularity at the grazing point and is given 

as 



X ^ ZDM{x) 



The correction S is given as 




if Hmin{x) > 
if Hmin{x) < 



(29) 



(5 = 2/^ 



{HxFi)xF 



v{x)^0{\xf) 



(30) 



where v{x) = Vi -\- V2 -\- vs with Vi^V2^vs G being each 
proportional to x and are given by 



vi 



V2 



-{- 



{HxF2)xF 
3 

{F2-Fi)xx, 



{F2-F,)xF- 



{Fi^xF2 - -Fi^xFi - -F2^xF2)xF} ^, 



^3 = -{F2-Fl)xF 



{HxF2)xX 

{HxF2)xF^ 



OH 

Hminix) = —X + 0{\xf). 



where the correction 7 is given as 

7 = Vi^{-H{x))i +V2zx{-H{x))^ 

V3AF2,:cX{-H{x)y^ +0{X^) 



(34) 



where 



Viz 



2 1 2 

jp^ jp^3/2 {'^i-^'^yxx ~ Fi^xx)F + F2^xFi^xF 

-\[Fl, + 2F2AF - „ I p (F2,x - Fi,.) 
F[-HxF2^xxF -\- HxF2^xFi^xF — -HxF2^xxF ]}, 



V2z 



2V2 



yHxFi^ 
2V2 



HxF2^xF^ HxFi^xF 



{F2,x — Fi^x)i 

{F2,x — Fi^x)F. 



(35) 



The subscript x denotes partial differentiation with re- 
spect to the variables x = {xi, X2, X3, X4} and the sub- 
script z denotes the projection along z, the normal vector 
to Eat. Further F denotes that the vector fields are con- 
tinuous at the discontinuity, that is Fi{x*) = F2{x^) = 
F{x^). Substituting Fi, F2, F, H{x) for the system de- 
fined by Eq. (25), the ZDM correction 5 as well as the 
PDM correction 7 for the grazing orbit at the transver- 
sal intersection with the discontinuity surface Si ^2 can 
be obtained. 

For the memristive MLC circuit, we define two scalar 
functions Hi{x) = (x — x*), = —1 and H2{x) = {x — 
X*), X* = +1, such that their zero sets define the two 
discontinuity boundaries J^i^2 and 112,3 • 

When the conditions for transversal grazings given by 
Eq. ( 26 ) are satisfied at each of the discontinuity bound- 



aries, then the ZDM correction 5 given by Eq. (30) can 
be derived by applying Eq. (25) in Eq. (31). 



Vi 



( —\{a\ — a2)xi 

|(ai - a2)(ai - 2(a2 + P))x^ + 



(31) 



V 



|(ai - a2)Pxi 




where 0(|xp is the correction term in the Taylor expan- 
sion of H{x). 

Let X be a point in ITat, the Poincare section defined 

by 

n^ = {xGro: ^|^Fi(x)=0}. (32) 

Then, sufficiently close to a regular grazing point satisfy- 
ing Eq.(26) of a second order non-smooth system defined 



V2 




(36) 



(37) 



above, the Poincare-section discontinuity map is given as 



X ^ PDM{x) = { ; n 

I X + 7 II H[x) < . 



'^3 




(38) 



Substituting these, the ZDM correction term S is given 



as 
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\V2{ai 



\V2{ai - a2)(ft)^/^(ai - 2(a2 + P))x3 + 2fsin{ujX4) 



\V2{ 



ai 



V 







(39) 

Similarly the PDM correction can be derived by sub- 
stituting Eq. (I25| & (l35|, in Eq. 



V2z 



2V2(ai-a2) 





2V2{ai-a2) 





(41) 



(42) 



Viz 



4(l+ai+^) I 4/ sin(a;a;4) 
3^3 ^3 3/^ 
2(ai-a2)(ai+/3) , 2(ai -a2)(ai+/3) . 

4(a^+a2+ai/3+/3^) 

6(ai— a2)(ai+/3)a;3+2( — 2+ai— 3a2 — 2/3)/ sin(u;g;4) 
3^ 

4/3(l+ai+/3-/3^) 4/3Vsin(u;a;4) 
3v^ ^ 3^ 





Making use of these the PDM correction term 7 is given 



as 



/ 

2(ai-a2)(ai+/3)xi3/=' 



V 



(43) 
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